Excitation spectrum of bosons in a finite one-dimensional circular waveguide via the 

Bet he ansatz 
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The exactly solvable Lieb-Liniger model of interacting bosons in one-dimension has attracted 
renewed interest as current experiments with ultra-cold atoms begin to probe this regime. Here 
we numerically solve the equations arising from the Bethe ansatz solution for the exact many-body 
wave function in a finite-size system of up to twenty particles for attractive interactions. We discuss 
the novel features of the solutions, and how they deviate from the well-known string solutions [H. 
B. Thacker, Rev. Mod. Phys. 53, 253 (1981)] at finite densities. We present excited state string 
solutions in the limit of strong interactions and discuss their physical interpretation, as well as the 
characteristics of the quantum phase transition that occurs as a function of interaction strength in 
the mean-field limit. Finally we compare our results to those of exact diagonalization of the many- 
body Hamiltonian in a truncated basis. We also present excited state solutions and the excitation 
spectrum for the repulsive ID Bose gas on a ring. 

PACS numbers: 03.75.-b, 03.75.Hh, 03.65.Ge 



I. INTRODUCTION 

Physics in low dimensional systems has long provided 
a rich source of fascinating and often unexpected phe- 
nomena. With the steady progress of experimental 
methods in ultra-cold gases, effective onc-dimensional 
systems are beginning to be realised in the labora- 
tory [U, i, i, il, i, 0, i, i, The integrability of 
certain many-body problems in one dimension pd] . [l^ 
provides an opportunity to reliably examine many-body 
quantum physics beyond mean-field theory. It is becom- 
ing clear that many condensed matter theories of phase 
transitions and collective excitations depend on the num- 
ber of degrees of freedom within the system [l^, . 

In this paper we focus on a system of identical bosons 
tightly confined in a ring trap such that the system can be 
considered to be purely one-dimensional (but with three- 
dimensional scattering). We are primarily interested in 
obtaining the excitation spectrum of the system. The 
motivation for this work was provided by the prediction 
of a quantum phase transition in the regime of attrac- 
tive interactions. As the interaction strength for a fixed 
number of particles becomes more negative, zero tem- 
perature quantum fluctuations eventually cause the gas 
to form a soliton-like localised state [13, [2l| . Quantum 
solitons in ID Bose gases were first predicted and ob- 
served with photons in optical fibres p^ . [l6t , and more 
recently soliton-like behaviour has been observed in sys- 
tems of massive particles with attractive interactions in 
quasi-lD [13, [Hi and 3D geometries [H]. 

The quantum phase transition in the ID toroidal geom- 
etry was initially identified using a mean-field approach 
by Kavoulakis [23| andKanamoto et. al. [2ll . l2^ . Further 



work by Kanamoto et al. has taken a quantum many- 
body approach using exact diagonalization by truncating 
the Hilbert space to three (or five) single-particle states 
for up to 200 atoms [H, [Ij] . They found a number of 
interesting physical features in the region of the phase 
transition, in particular evidence of symmetry-breaking 
as the difference in energy between the ground and first 
excited stated tended to zero and scaled as iV~-^/^. 

Whilst the results obtained by Kanamoto et al. [H, 
remain qualitatively correct at low temperature and weak 
interactions, the Bethe ansatz provides us with an exact 
solution for all excited states of the system at all values 
of interaction strength. The ground state quasi-momenta 
and energy for the onc-dimensional Bose gas has previ- 
ously been calculated using the Bethe ansatz by Sak- 
mann et al. (25j . However, we emphasize that for exper- 
imentally realistic temperatures for this system excited 
states will be involved in both the static and dynamic 
properties of the system. In this paper we extend the 
work of Sakman et al. [25| to calculate the low-lying ex- 
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citation spectrum of the finite one-dimensional Bose gas 
as a function of interaction strength for up to A'^ = 20 par- 
ticles. This is a non-trivial result due to the complicated 
behaviour of the quasi-momenta in the complex plane in 
the attractive case (see section UlI Ap . To the best of our 
knowledge this work is the first example of calculations 
made on excited states of this system (with more than 
3 particles) without any truncation of the Hilbert space 
and without restricting oneself to the limiting string so- 
lutions. 

The string solutions arise in the case where the inter- 
atomic interactions are sufficiently attractive, or alter- 
natively in the zero density limit. The problem reduces 
to one solved by McGuire [l2| and elaborated upon in 
references [26l [27j . The point of interest in our numeri- 
cal solutions is the deviations of the quasi-momenta from 
these previously known solutions. Our present work ex- 
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tends our understanding of quantum solitons from the 
large boson numbers in optical fibres towards the much 
smaller numbers possible in atomic systems. 

Related work has been performed by Oelkers and 
Links [2^, who concentrate on a toroidal lattice gov- 
erned by the Bose-Hubbard Hamiltonian with periodic 
boundary conditions. Physically one expects the results 
of their calculations to tend toward those of the contin- 
uum as the number of lattice sites increases. The Bose- 
Hubbard Hamiltonian is, however, non-integrable and it 
has been shown by Seel et al. [2^ that the spin-1/2 XXZ 
Heisenberg chain (which is integrable) maps onto the ID 
Bose gas model in the continuum limit. Other examples 
of lattice Hamiltonians that have the ID Bose gas in the 
continuum limit can be found in [s^, [3l| . 

The structure of this paper is as follows. In section HIl 
we review the work of Licb and Liniger pd| in order to 
correctly pose the problem. In section IIII Al we discuss 
our method of solution and present the results, whilst also 
discussing its connection to earlier work by McGuire . 
In section UlI BI we present the excitation spectrum of the 
system obtained from our calculations, and in section 
IIII CI we show a comparison of our work to that of the 
truncated Hilbert space approach of reference [2^, [2^ . 
Finally we calculate the excitation spectrum for the re- 
pulsive case in section ITVl before concluding in sectionlVl 



II. THE MODEL 

We are considering a system of N bosons, each with a 
mass m, on a ring of radius i?, in the limiting case of the 
radial confinement freezing out all the transverse degrees 
of freedom of the atoms. The interactions are assumed 
to be short range, contact interactions, modelled using 
the delta function. With these assumptions the system 
will be described by the Hamiltonian (in first quantised 
formalism) 
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where length is measured in units of R and energy is mea- 
sured in units of /2raR^ {Oi is the angular coordinate 
of the ith particle). The parameter c is related to the 
s-wave scattering length by Eq. and quantifies the 
strength of the two body interactions, c < describes an 
attractive gas while c > describes a repulsive gas. 

Since our spatial universe is of course three dimensional 
it is worthwhile to note the physical limits of applicabil- 
ity of this model. Assuming harmonic confinement, the 
length scale of the dimensions tranverse to the ring are 
xq = \J hj (miOx) and yo = \/hJ(mxJy). Another impor- 
tant length scale is the s-wave scattering length, a, that 
represents the strength of two body interactions. The 
healing length of the gas ^ — l/y/Snna, where n is the 
three dimensional density, reveals how sensitive the gas 
is to irregular features in the potential. For the system 



in which we are interested, it is possible to define a local 
set of cartesian coordinates such that zq 3> C ^ xo,yo, 
where the z axis runs tangential to ring. Essentially this 
amounts to a torus of radius R and cross-sectional area 
7rr^ satisfying the condition R ^ r. In this way the en- 
ergy level spacings for the tranverse (x and y) degrees of 
freedom supersede all other energy scales within the sys- 
tem and the system is essentially one dimensional. There- 
fore we can assume a ground harmonic oscillator state in 
the transverse dimensions, and integrate them out to ob- 
tain an effective one dimensional interaction 13211 



(2) 
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where to is the mass of a single particle. Some interes ting 
work has recently been done by Parola et. al. [13, l33 | 
on the quasi-one dimensional limit which indicated a 
transverse collapse (in addition to the angular collapse) 
at a certain critical interaction strength. In this work 
they took into closer consideration the effects of the ex- 
cluded two dimensions. Current experimental realisa- 
tions of such ring t rap s 13511 are still far from the quasi- ID 
regime [H, 113, S HI, 

Given the Hamiltonian ([1]), Schrodingcr's equation is 



,0Ar) = K^„(ei,...,ejv), 



(3) 



where n is the label for different eigenstates. Equation ([3]) 
can be solved exactly via the Bethe ansatz Due to 
the symmetry under permutation of particle coordinates, 
the region over which we need to integrate Eq. ([3]) can be 
restricted to < 0i < . . . < ^at < 2tt. Lieb and Liniger 
noted that in this region the delta function vanishes ev- 
erywhere except along the boundaries where 9i = 6i^i, 
hence the problem could be cast in a different manner 
by writing the interactions as a boundary condition on 
the restricted reg ion rather than an explicit term in the 
Hamiltonian The new problem was to find the ap- 
propriate solution to 

— J2jLi ^ai^ ~ E„il'n in the restricted region, (4) 



(5) 



where Eq. ([5]) is the previously mentioned boundary con- 
dition arising from the interaction. Furthermore the 
boundary conditions pertaining to the ring geometry of 
the system, i.e. ipni- ■ ■ , 6*1, • ■ ■) = ■>Pn{- . . , 6*^ -I- 27r, . . .), 
must also be recast onto the new region, vis 



27r), (6a) 



do 
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The Bethe ansatz is employed as a means of integrat- 
ing Eq. (j4]). The (unnormalized) wave function in the 
restricted region is given the form 



V'n(^l,---,^JV) 



{Q} 



(7) 
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where Q is some permutation of the integers 1, 2, . . . , iV, 
Q{j) is the jth element in Q and the sum runs over all 
the A^! distinct permutations of Q. Equation ([7|) is one 
form of the Bethe ansatz and was originally used as a 
means of integrating spin chain Hamiltonians [4l|. The 
coefRcients of each permutation, A^q^ , can be determined 
(up to some arbitrary phase) by the interaction boundary 
condition Eq. ([5]). We choose the convention ^^2' n = 
1 and for any other permutation Q, Aq''^ is found by 
constructing Q out of 12 ... by transposing adjacent 
elements. Then Aq'^ is a product of the terms 



J") 
«Q(b) 



J") 



Q(b) 
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where Q{a) and Q{b) are the adjacent elements being 
transposed, with Q{a) to the left of Q{b) before the trans- 
position. 

The scalar quantities /cj"^ , • ■ • , fcj^''' , known as the 
quasi-momenta (and sometimes the rapidities), are de- 
termined by the ring geometry of the system, i.e., by the 
boundary conditions in Eqs. These quasi-momenta 
furnish important quantities for the system such as the 
energy eigenvalues and the total momentum, however it 
should be stressed that (other than the case of no in- 
teractions c = 0) the individual quasi-momenta do not 
have a direct relationship to the single particle momen- 
tum states. The latter must be obtained from the single 
particle reduced density operator, which can be obtained 
from the many body wave function. 

Substitution of Eq. ([7]) into the boundary conditions 
((6a|) and ((6b)) . yield the following set of TV simultaneous 
equations to determine the quasi-momenta as functions 
of c 
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N 

n 
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-, j^l,...,N. (8) 



The purpose of this work is to explore solutions of these 
equations to reveal key features of the system such as the 
excitation spectrum of the g function of the inter- 

action strength. Analytical solutions of Eq. ([5| have pre- 
viously been found for up to three particles |42l | , however 
to our knowledge the only numerical calculations have 
been for the ground state energy Once Eqs. ([8]) have 
been solved it is true that, in principle at least, one has 
obtained all neccessary information to obtain the wave 
function for the system. However in practice, even for a 
modest number of particles such as ten, the wave func- 
tion itself will involve a summation over 10! ~ 4 x 10^ 
terms, a somewhat cumbersome object. However, for 
certain physical quantities, simple analytical expressions 
involving only the quasi-momenta can easily be derived, 
for instance the energy of the nth eigenstate is given by 



E„ 



N 



(9) 



and likewise the total momentum by 

N 



(n) 



(10) 



We will utilise equation ([9|) in section UlIBI to obtain the 
excitation spectrum. 



III. THE ATTRACTIVE GAS: c < 

The behaviour of the Bose gas with attractive inter- 
actions has received far less scrutiny than its repulsive 
counterpart. Originally Lieb and Liniger did not con- 
sider this regime to be of physical relevance [ll|, how- 
ever McGuire showed [l2| that once an A^-particle bound 
state had formed, the ground state energy of the system 
scales as 



Eq cx -c^N{N^ - 1). 



(11) 



This posed problems for the solution in the thermody- 
namic limit where N —i- 00, R 00 and N/R = con- 
stant. In contrast for the repulsive case it is found 
Eq oc N (X R, and hence an expression for the energy 
density could be found [H, 0, 0, However E/R 
diverges in the attractive regime. It is perhaps worth 
noting that the limit of an extremely dilute gas is free 
from these divergences, provided one took the limit as 
/ R = constant. Another possibility is to consider the 
limit R^ 00 while N remains constant [H, [13, [13, [ill . 

As previously mentioned the quantum phase transi- 
tion to a soliton-like state provides added interest for 
the attractive gas. This transition spontaneously breaks 
the translational symmetry of the system even at zero 
temperature. We discuss this point further in section 
IIII Al and see how the eigenstates obtained from the 
Bethe ansatz maintain the translational symmetry of 
the Hamiltonian by forming a continuous superposition 
of localised states around the ring. We have included 
in appendix |B] a simple derivation (following the work 
of Kavoulakis of the critical interaction strength 

(Co = —tt/2N) at which the phase transition is predicted 
to occur. 



A. The Quasi- Momenta 

To find the quasi-momenta {ki} that fully characterise 
the exact solution, we must numerically solve Eqs. ([5]). 
For N particles this gives N simultaneous nonlinear equa- 
tions. The nature of most root-finding algorithms is such 
that it requires a reasonable initial guess for the location 
of the roots. We have developed a suitable procedure to 
solve these equations for a range of interaction strengths 
c starting from the known solutions for the ground and 
excited states at c = 0, which arc given by all ki being 
an integer. 



4 



(a) 



Id 0.5 




Id 0.5 




0.5 



V 0.5 




(d) 



Id 0.5 









Id 0-5 




Id 0.5 




Id 0.5 




(g) 



(h) 



0.5 






Id 0.5 



Id 0.5 



FIG. 1: (colour online) The behaviour of the quasi-momenta for A'^ = 6 particles as a function of the interaction strength 
c < 0, for the first 8 states, ordered via their energy at c = 0. On each graph there are six lines, the color of which indicates 
how many quasi-momenta are overlapping. Blue indicates one root, green indicates two roots, magenta indicates three roots, 
cyan indicates 4 roots and red indicates six roots. The single particle momentum states at c = are, (a) {0,0,0,0,0,0}, 
(b) {1, 0, 0, 0, 0, 0}, (c) {1, 1, 0, 0, 0, 0}, (d) {1,-1, 0, 0, 0, 0}, (e) {1, 1, 1, 0, 0, 0}, (f) {1, 1, -1, 0, 0, 0}, (g) {1, 1, 1, 1, 0, 0} and (h) 
{1, 1, 1, —1,0,0}. The states with a non-zero total angular momentum (that is states (b),(c),(e),(f),(g) and (h)) have all been 
chosen to have positive total angular momentum and naturally they have a degenerate counterpart with a negative total angular 
momentum which can be found easily via the mapping ki —ki for all i. 



To find the solution for non-zero c, we first choose a 
value for c close to zero, and use the ideal gas solution 
for the state we are interested in as the initial guess for 
the root-finding algorithm. We make use of the built-in 
root-finder f solve in the software package MATLAB on 
a standard desktop PC. By choosing a small enough c 
the initial guess for the quasi-momenta are close to the 
actual solution, and the root-finding algorithm converges 
relatively rapidly [55| . 

Once the first solution close to c = is found, c is de- 
creased in small increments Ac and an initial guess for 
the quasi-momenta for the new value of c is based on 
a smooth extrapolation from the previous values. Typi- 
cally the extrapolation is based only on the two previous 
values of c and is thus linear. Difficulties arise in this 
method when a particular quasi-momenta at an interac- 



tion strength c differs greatly from that at c -|- Ac. This 
problem arises particularly when c is close to the critical 
interaction strength cq and we discuss how we deal with 
the problem later. 

The actual size of Ac will be dependent on how large 
N is. The scaling behaviour of the critical interaction 
strength goes as Co (x 1/N and wc find that the efficiency 
of the algorithm is improved by using Ac as small as 
10~^/N. It is possible to use a larger value of Ac but 
wc have found that smaller step sizes give better initial 
predictions for the roots for the next step and results in 
faster convergence of the root-finding algorithm. 

We have not found it necessary to make use of high pre- 
cision arithmetic used by Sakmann et al. [25t in finding 
our solutions. In their work they make use of a particular 
transformation [see Eq. (33) in Ref. [2^] popularised in 
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the seminal work of Lieb and Liniger [ll| . The transfor- 
mation simpUfies the numerics in the repulsive case by 
spreading the transformed quasi-momenta into numbers 
which can be easily distinguished by machine precision. 
The same is not true in the case of attractive interac- 
tions. As pointed out by Sakmann et ai, working with 
the transformed quasi-momenta in the attractive case re- 
quires numerical precision of approximately 10~*^. In 
this work we solve Eqs. ([8]) directly. We set the toler- 
ance of the f solve algorithm to iterate until Eqs. ^ arc 
solved such that the left hand side equals the right hand 
side to an absolute accuracy of 10~^°. Once the value of 
|c| has been increase to be point that the quasi-momenta 
are within 10~^° of the string solution the algorithm ter- 
minates. 

We present the results of these calculations in two dif- 
ferent ways. Figures [U [21 and [3] show the continuous 
evolution of the real and imaginary parts of the quasi- 
momenta, for a number of excited states as c becomes 
more negative. Figure [1] shows eight excited states for 

= 6, whereas Figs. [2] and [3] show two representative ex- 
cited states for = 20. Figures [H and [5] show snapshots 
of the distribution of the quasimomenta in the complex 
plane for A^ = 20. Our calculations elucidate the devia- 
tions from the string solutions at finite density (49l | . 

We see from Fig. [TJb) and (e) that as the real parts 
of ki and kj become equal the imaginary parts bifurcate, 
and vice-versa. These splittings of the quasi-momenta 
provide additional challenges for the root-finding algo- 
rithm. Because the splitting often occurs sharply, the 
initial guess for the quasi-momenta found from the pre- 
vious value of c can be sufficiently inaccurate that the 
algorithm does not converge. Thus if one finds that the 
root-finding algorithm is not converging at a particular 
value of c then this could be an indication that two equal 
quasi-momenta are beginning to split. When this occurs 
there is always a degeneracy between the quasi-momenta 
as to which one goes up and which one goes down (this 
degeneracy is of no physical consequence, it is merely a 
mathematical hurdle) . The numerical fluctuations of the 
values of the quasi-momenta is usually somewhere close 
to machine precision, w 10~^^, and are difficult to detect. 
In order to observe the splitting we found it was necessary 
to manually alter the initial guess for the quasi-momenta 
at these points, forcing one to go up and one to go down 
(again see Fig.[T](b) and (c)). 

There is no fundamental limit to the excited states we 
can reach with this procedure, and similar graphs can be 
produced for up to approximately A^ = 20 particles (see 
Fig. [2] and [3] as examples of the A^ = 20 particle gas) . 
In principle this work could be extended to larger num- 
bers of particles, and in fact we have found the ground 
state for up to A^ w 50. However, the bifurcations in the 
excited states require some time to locate and we have 
had to deal with them manually. This means that the 
algorithm begins to take a considerable amount of time 
for A^ significantly larger than twenty. 

Roots of Eq. ([S]) for c < have a far more complicated 
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FIG. 2: Quasimomenta of the first excited state for A = 20 
particles as as a function of the interaction strengnth c. (a) 
and (b) show the real and imaginary parts respectively of the 
quasi-momenta. The blue solid lines are the roots calculated 
using the Bethe ansatz, whereas the red dashed lines show the 
string solutions given by Eqs. (|13p which give the asymptotic 
behaviour of the quasi-momenta in the attractive limit. The 
two solutions are in good agreement above the critical interac- 
tion strength, validating the essential pro perties of the bound 
states described in references p^. [27l,l49l. [sol. IsH. The inset 
shows the region where the quasi-momenta bifurcate causing 
numerical difficulties. 



behaviour than in the more commonly studied repulsive 
regime. Apart from the obvious reason that the quasi- 
momenta can be c omp lex (as opposed to strictly real for 
the repulsive gas [4^), there also exists the change in 
behaviour at the critical interaction strength Cq. For 
some time it has been suspected that the roots for large 
enough attraction would take values on strings in the 
complex plane which correspond to some bound state of 
the system [i^, HO, [Ml ■ That is, the wave function would 
represent some kind of localised state [l^, such as 

V'6ow(0i,...,0iv) =exp I -i|c|^|0,-0j| , (12) 
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FIG. 3: Quasimomenta of one of the second excited states for 
N = 20 particles as as a function of the interaction strengnth 
c. (a) and (b) show the real and imaginary parts respectively 
of the quasi-momenta. The blue solid lines are the roots cal- 
culated using the Bethe ansatz, whereas the red dashed lines 
show the string solutions given by Eqs. (|16|l with M = 1 The 
inset shows where the quasi-momenta bifurcate for this case. 



and the corresponding quasi-momenta of this state would 
be 

K 1 , 
ki = - + -iN-l)tc, 



(13) 



where K — ki is the total momentum of the state. 
It is worthwhile to note that although we refer to the 
wave function in Eq. (jl2p as a localised state it is not 
localised to any specific point. Rather it has a localised 
pair correlation [47[ 



(*t(6))#(6i))(^t(6)')#(6i')) ' 



where ^{9) {^'^{9)) are bosonic field operators which an- 
nihilate (create) a particle at position 9. Thus measure- 
ments of the density of an ensemble will always yield a 
localised projection, whilst the system is not localised 
prior to the measurement. One can think of the system 
prior to the measurement as being in a (macroscopic) 
superposition of localisation at every point on the ring. 
This form for the quasi-momenta given by Eqs. p3p is 
commonly referred to in the literature as a string solu- 
tion [i^, [H^, HH . It is often assumed that multiple strings 
may exist for one system ([i^ and references therein), 
with each string corresponding to a soliton of different 
momentum. This choice of quasi-momenta also seems 
intuitively reasonable and gives the result, mandated by 
McGuire for the total energy of the bound state 



I NjN^ 1) 
^''^N^ 12 



(15) 



(14) 



The degree to which the system is localised by these 
string solutions can be quantified by the absolute val- 
ues of the imaginary parts of the quasi-momenta. Thus, 
in Eqs. as |c| increases, so does the degree of locali- 
sation (as one would naively expect). However, it is easy 
to show that in the ideal limit, c — > 0, the quasi-momenta 
are exactly the single particle momentum states. Thus if 
one accepts the quasi-momenta given by Eqs. then 
one is left with chronic discontinuities in the behaviour 
of the quasi-momenta as a function of |c|. Our results 
give direct evidence that the two limits are bridged via 
Eqs. (IHl). We also find other families of excited states, 
with multiple bound states forming in the limit of strong 
attraction. These show agreement with the truncated 
diagonalization results discussed later in section IIII CI 
These families of excited states can be interpreted as the 
formation of multiple bound states. It is very interesting 
to observe this intermediate behaviour between the ideal 
gas and strongly attractive gas since this is where soliton 
formation occurs at zero temperature. 

In Fig. [T] we see the behaviour of the quasi-momenta 
for a system of iV = 6 particles in the eight lowest energy 
states ordered via their energies at c = 0. As |c| becomes 
larger, energy level crossings are observed as in Fig. [S] 
The quasi-momenta begin at c = at their respective 
single particle momentum states, and as c becomes more 
negative the quasi-momenta tend toward a string solution 
as predicted. The initial behaviour of the imaginary part 
of the quasi-momenta goes like ^/\c\ for |c| < Co, but 
for |c| > Co when the system is in a localised state the 
low-lying eigenstates have quasi-momenta corresponding 
to the string solutions Eqs. (fT5|) . 

This is verified again for the = 20 particle gas in 
Fig. [2] where the solid blue lines show the quasi-momenta 
as a function of |c| and the red dashed lines show the 
string solution Eqs. p^ . and both are in agreement in 
the strongly attractive limit. However, we also find ex- 
cited states that do not agree with Eqs. (fTS]) . but instead 
form multiple bound states. This results in a decrease 
in the degree of localisation. In these cases we find the 



7 



(a) 



(b) 



(c) 



c- 



c=c 



c=1.9c 



Re(k,) 



Re{k.) 



Re(k.) 



FIG. 4: Quasi-momenta distribution for A'^ = 20 particles in the complex plane for different values of interaction strength c for 
the first excited state as in Fig. (2] (a) |c| = 10^^. The quasi-momenta are close to the single particle momenta for the ideal 
gas. (b) |c| = Co ~ 0.079. Near the mean-field critical point there exists a two-particle bound state, a three-particle bound 
state and a fifteen-particle bound state, (c) |c| — 1.9. The interaction strength is well past the critical point and the three 
particle bound state has collapsed in with the fifteen particle bound state, however there is still a two particle bound state. 






FIG. 5: Quasi-momenta distribution for A'' = 20 particles in the complex plane for different values of interaction strength c for 
the first excited state as in Fig. [31 (a) |c| = 10"^ The quasi-momenta are close to the single particle momenta for the ideal 
gas. (b) |c| = Co ~ 0.079. Near the mean-field critical point, there exists a free particle, a five-particle bound state and a 
fourteen-particle bound state. (c)|c| = 1.9. The interaction strength is well past the critical point and the five-particle bound 
state has combined with the fourteen-particle bound state leaving a single free particle. 



behaviour of the quasi-momenta at large |c| to be given 

by 



— + -(M~l)ic 



K 1 

M ~ 2 
K' 



N -M 2 



(M - l)ic 
+ -{N ~M - \)ic 



K' 



N ~M 



M~l)ic. (16) 



We interpret this solution physically as an M particle 
bound state with momentum K and an iV — M particle 
bound state with momentum K' . 

Evidence of these states can be seen in Fig. [3] where the 
solid blue lines show the quasi-momenta as a function of 
|c|, whereas the red dashed lines show the solutions given 



by Eqs. ^ with M = l,K = -1/20 and K' = 21/380. 
In Fig. [Ijd) and (f) the same kind of behaviour is seen 
for the iV 6 particle gas. The behaviour is perhaps 
clearer in Figs. (DJc) and OJc) where one can see the two 
distinct strings in the complex plane. One can derive a 
simple expression for the energy of these states (similar 
to Eq. (HID), 



E 



(1) 

KK' 



M{N - M) 
N{N^ -iNM + MP - I) 
12 



(17) 



This trend continues such that three or more bound 
states of atoms occur within the system. The evidence 
for this is most easily seen by the grouping of the quasi- 
momenta in the complex plane [Figs.Hfc) and[5{c)]. 

Thus for this system once c < Co and localisation has 
occurred we observe several different families of solutions. 
The first is the ground state of the system, corresponding 
to a single, stationary A^-particle soliton. Then there are 
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FIG. 6: (color online) The excitation spectrum of a A'^ = 
6 particle system obtained via the Bethe ansatz. The solid 
red lines show the formation of a single bound state. The 
green dashed lines show the formation of two bound states. 
The blue dot-dash lines show the formation of three separate 
bound states. The qualitative change in the behaviour of the 
eigenstates at the phase transition point becomes more visible 
than for the two particle case (see Fig. [Sjl . 

the elementary excitations of this state whereby the soli- 
ton has some (integer valued) total momentum about the 
ring. Furthermore there exist the higher order excitations 
in which multiple solitons form around the ring. Equally 
these multiple solitons can have elementary excitations 
of their own, corresponding to some integer valued total 
momentum. 



B. The Excitation Spectrum 

We can now calculate the excited states of the system 
via the Bethe ansatz, using Eq. ([9]) and the numerically 
determined roots of Eqs. ((8|) (the quasi-momenta). We 
show a comparison between the truncated diagonaliza- 
tion approach of Kanamoto et al. [23l. [23| in section UlI CI 
(sec Figs. [8] and [9]). Wc plot the excitation spectrum for 
iV = 6 particles in Fig. [HI and = 20 particles in Fig.[71 
The point at which mean-field theory predicts the quan- 
tum phase transition, Co (see appendix iBjl . is indicated 
on all figures with a vertical dot-dashed line. We observe 
the sharpening of the cross-over regime as N increases as 
expected for a phase transition in the limit of large N. 

Once the size of the system gets up to around = 20, 
the difference in the behaviour of the states for |c| < Co 
and |c| > Co is clear. For |c| > Co transition (once the 
bound state has formed), the different families of solu- 
tions to the Bethe Eqs. ([8]) (discussed in the previous 
section) become distinguishable through their separation 
from the ground state. The solid red lines in Figs. [6] and 
[7] show the behaviour of the single bound states found by 
McGuire and others [H, 113, Eo, [5l| given by Eqs. (fTS]) . 
the green dashed lines show the family of solutions given 



FIG. 7: The excitation spectrum of a A'^ = 20 particle system 
obtained via the Bethe ansatz. The solid red lines show the 
formation of a single bound state. The green dashed lines 
show the formation of two bound states. The blue dot-dash 
lines show the formation of three separate bound states. 



by Eqs. (|I6p which are made up of two separate bound 
states. Finally the blue dot-dashed lines show the exis- 
tence of three separate bound states within the gas. 



C. Comparison to the truncated Hilbert space 
diagonalization 

The excitation spectrum of the attractive gas has been 
numerically studied via truncation of the Hilbert space 
in previous work [23l . [23 |. The key element in the ap- 
proach is to exclude all single particle momentum states 
beyond a particular threshold. The results of these cal- 
culations are of course exact at zero interactions, but de- 
viate from the exact solution as the interaction strength, 
|c|, is increased. The degree of deviation can of course be 
decreased by increasing the threshold momentum state, 
however the dimension of the Hilbert space will eventu- 
ally become unmanageably large in doing so. An example 
of this process is shown in Fig. [5] for a two particle sys- 
tem where it is possible to make the cut-off relatively 
high without having the Hilbert space dimension expand 
beyond a practical number. The results indicate a good 
qualitative description from a very small number of mo- 
mentum states [231 . For example if A'^ = 200, a momen- 
tum cutoff of ±1 is reasonable up to and even beyond 
(c « 2co) the critical point. This was supported by cal- 
culations made by Kavoulakis [2^ who used a mean field 
approximation with a suitably chosen wave function to 
examine the ground state and low-lying states near the 
phase transition. 

The truncated diagonalization procedure we use fol- 
lows that of Kanamoto et al. [ISl- The first quantised 
Hamiltonian in Eq. ([T|) is rewritten in second quantised 
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FIG. 8: Comparison of the excitation spectrum found via 
the Bethe ansatz Eqs. ([8]) (red dashed lines) to that found 
via diagonalization of the truncated Hamiltonian (solid blue 
lines) for a, N — 2 particle system. The vertical dash-dotted 
line marks the mean field crititical interaction strength Co 
derived in appendix [B] (a) Single particle momentum state 
cut-off of ±1. (b) Single particle momentum state cut-off of 
±2, and so on, up to (f) which has a single particle momentum 
cut-off of ±6. 



formalism as, 

n27r 

H 



de 







(18) 



where ^(0) is the bosonic field operator that annihilates a 
boson at coordinate 9 and obeys the usual commutation 
rules and the periodic boundary conditions of the ring 
geometry. This field operator is then approximated by 
its truncated expansion in terms of single particle states 



(19) 



- — ko 



where Lpj (9) is the single particle eigenstate with angular 
momentum j (see appendix VK\ . Cj annihilates a boson 



FIG. 9: Comparison of the excitation spectrum found via 
the Bethe ansatz Eqs. ([SJ (red dashed lines) to that found 
via diagonalization of the truncated Hamiltonian (solid blue 
lines) for a A'^ = 20 particle system, (a) Truncation at m = 3 
(single particle momentum states 0,±1). (b)Truncation at 
m = 5 (single particle momentum states 0, ±1, ±2). 



with angular momentum j and fco is the single particle 
momentum state cut-off. In this manner, for finite N and 
finite fco 5 -ff is a finite matrix. Furthermore, H (when 
written in the appropriately ordered basis set) will be 
block diagonal, as it can be divided into total momentum 
subspaccs which can be individually diagonalised. If we 
use m to denote the total number of single particle mo- 
mentum states used in the expansion (i.e. m = 2fco + 1) 
then the total dimension of the Hilbert space will be the 
binomial coefficient (^+™~^). This Hilbert space will 
split up into N{m — 1) + 1 (total momentum) subspaces. 
The largest subspace will be the zero momentum sub- 
space. For example if we were to consider a gas of = 20 
particles with the momentum cut-off at ±2 (m = 5), then 
the total dimension of the truncated Hilbert space would 
be 10626. However this splits up into 81 different total 
momentum subspaces, the largest of which has dimension 
318, which can be exactly diagonalised on a standard 
desktop PC. Table |T] gives an indication of the scale of 
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81 


318 
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121 


5444 


200 
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20301 


401 


101 




5 


70058751 


801 


230673 



TABLE I: The size of the Hilbert space and subspaces for 
selected particle numbers and momentum cut-offs. This gives 
an idea of the size of computation required for the truncated 
diagonalisation procedure. 



computation required for the truncated diagonalization 
procedure. 

In Fig. [9] we show a comparison between excitation 
spectrums obtained from the truncated diagonahzation 
procedure and the Bethe ansatz for the TV = 20 particle 
gas. We observe that for m = 3 (single particle momen- 
tum cut-off at ±1) the truncated diagonalization over- 
estimates the energy difference between ground and ex- 
cited states in the cross-over regime, but undcr-cstimatcs 
the difference in the regime where the bound states have 
formed. For to = 5 (single particle momentum cut-off 
at ±2) the agreement between the two methods in the 
cross-over regime is vastly improved, however the trun- 
cated diagonalization still underestimates the difference 
between ground and excited states in the bound state 
regime. 



IV. THE REPULSIVE GAS: c> 

The root-finding algorithm can easily be applied to the 
case of repulsive attractions. A proof exists showing that 
all roots of Eqs. ^ lie on the real axis [i^. This means 
that we avoid the problems associated with the bifur- 
cations which appeared in the attractive case and the 
situation is therefore significantly easier. In this case the 
limiting behaviour as c ^ oo of the quasi-momenta is 
given by the siiigle particle momentum states of an ideal 
Fermi gas [ill, 113 ■ This limit is known as the Tonks- 
Girardeau gas and has been experimentally observed ^ . 

In this work wc use our numerical root-finding algo- 
rithm to bridge from the ideal gas to the Tonks-Girardeau 
gas. We do so for both the ground state and the low-lying 
excited states for finite numbers of particles. The ground 
state has appeared previously in the work of Sakmann et 
al. psj . The thermodynamic limit has been addressed 
in the homogeneous case [i^ 13, and the trapped 
case [15 (and references therein). 



A. The Quasi- Momenta 

The root-finding algorithm is essentially unchanged 
from the case of attractive interactions. It is worth men- 
tioning however that, since there is no phase transition in 
the repulsive gas the step size Ac can be increased with- 
out sacrificing a reasonable initial guess. In our calcula- 
tions we use Ac = 10^"* for < c < 1, then Ac = 10"'^ 
for 1 < c < 6, then Ac = 10"^ for 6 < c < 20, then 
Ac = 10^1 for 20 < c < 100 and finally Ac = 1 for 
c > 100. 

In Fig. [To] we see how the quasi-momenta behave as a 
function of the interaction strength c. The single particle 
excitations of the ideal bosons at c = result in single 
particle excitations of the ideal fermions at c ~ oo. The 
concept of particle-hole excitations for the strongly repul- 
sive gas (see references [13, S^) is clear in Fig. [TUl We 
see low-lying excited state quasi-momenta are similar to 
that for the ground state, but with a small number of 
quasi-momenta excited, leaving behind holes. 

It is possible to push the algorithm for the repulsive 
case to much higher atom numbers than with the at- 
tractive case. This is due to the absence of the exotic 
features in the quasi-momenta spectrum (sec for exam- 
ple the splitting occurring in the insets of Figs and [3]) . 
In Fig. [11] we show the ground state quasi-momenta for 
a repulsive gas with N = 50. 

B. The Excitation Spectrum 

As for the attractive gas, the excitation spectrum can 
be found using Eq. In Fig. [T^] we show the ener- 

gies of the low lying excited states for c > 0. The red 
dashed lines show the result obtained from truncating 
the Hilbert space down to five single particle states and 
diagonalising the Hamiltonian (see section IIII Cp . The 
truncation procedure is accurate for small enough values 
of c (see the inset of Fig. [T^] where c runs from zero to 
one) , but fails as c increases beyond one. The dotted line 
shows the ground state energy of a TV = 20 particle Fermi 
gas. 

V. CONCLUSIONS 

We have studied the behaviour of the roots of the Bethe 
ansatz Eqs. ([H]) for the ID Bose gas on a ring as a func- 
tion of the interaction strength for both the attractive 
case (up to = 20 particles) and the repulsive case. 
We used these roots to determine the exact many-body 
energy spectrum of the one-dimensional Bose gas with S- 
function interactions in a ID ring trap for small numbers 
of particles. We compared our results to those obtained 
via the approximate method of truncated diagonaliza- 
tion, to which wc found good qualitative agreement and 
reasonable quantitative agreement for small enough in- 
teraction strength. Wc have quantitatively addressed the 
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FIG. 10: The quasi-momenta as a function of interaction strength for c > 0. (a) Ground state of the ideal Bose gas. (b) First 
excited state of the ideal Bose gas (degenerate with the state kj —kj)- (c), (d) Degenerate second excited states of the 
ideal Bose gas (also degenerate with the state kj —kj). (e), (f) degenerate third excited states of the ideal Bose gas (also 
degenerate with the state kj — > —kj). (g), (h) , (i) Degenerate fourth excited states of the ideal Bose gas (also degenerate with 
the state kj — > —kj). The insets on (a)-(i) show the behaviour close to c = 0. 



issue of when it is reasonable to assume a bound state 
solution to Eqs. ([8]), and have quantified the deviations 
from the string solutions that occur for a finite density 
gas, thus establishing the point at which the string solu- 
tions become valid. Furthermore wc found evidence for 
the existence of multiple, independent solitons existing 
on the ring. These excitations had been hypothesised in 
previous work, however to the best of our knowledge this 
is the first direct evidence of their existence. We have de- 
scribed the analytical expressions for the quasi-momenta 
of these new families of string solutions in Eqs. (|16p ) and 
an expression for the energies of these states. It may 
be interesting in the future to make use of the solutions 
found in this paper to calculate correlation functions for 
the one-dimensional Bose gas. 
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APPENDIX A: SINGLE PARTICLE IN A RING 

The single particle Hamiltonian is 



ifie) ^ Eipie), 



(Al) 
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c 

FIG. 11: The quasi-momenta as a function of interaction 
strength for an A'^ = 50 particle repulsive gas. 



with corresponding cigcnstates are 

^,(9) = -^e^^^ (A3) 

APPENDIX B: THE MEAN-FIELD CRITICAL 
POINT 

We include a derivation of the mean-field criti- 
cal point of localisation using the same approach as 
Kavoulakis [2^. Considering the Hamiltonian we 
define a mean field tp{0) and approximate the system 
as having the many body wave function Jli'/'(^i)- We 
next expand "0(0) in terms of single particle states given 
by Eqs. (jA3[) . however we truncate all momentum states 
greater than one, 

ijiO) = a-iip-i{0) + ao(po{9) + aiipi{e). (Bl) 



100 



LU 

I 

c 

LU 




100 



We then assume |ao| ^ |tt-iUcn| and from the symme- 
try |a_i| = |q!i|, and finally the normalisation condition 
|a_ip -I- |aoP + = 1. The energy per particle e is 
then found to be 



2\a,\ 



[l«o| 



-h4|aonaiP + 4|a-i 
+2{al)*a-iai] , 



2| 



4|ao 



|ai| 



(B2) 



e IS mm- 



where 7 = Nc/tt. Setting aj = \aj\e'^'^^ , then 
imised by choosing (j)^i + (f)i — 2(f)o = 0, the normalisation 
condition further simplifies Eq. (|B2p to 



FIG. 12: The excitation spectrum of a AT = 20 particle gas 
with repulsive interactions c > (solid black lines. The red 
dashed lines were obtained by truncating the Hilbert space 
down to five single particle momentum states. The inset 
shows the region of < c < 1 where the truncated diagonal- 
ization is accurate. The dotted lines shows the ground state 
energy of a free Fermi gas, and the so-called Tonks Girardeau 
ground state asymptotes to this value. 



I =2|ai|2(l + 27)- 77|ai|^ 



(B3) 



From Eq. (|B3p it is straightforward to see that the ap- 
proach one must take in order to minimise e will depend 
on whether — 1/2 < 7 < (in which case the gas is delo- 
caliscd) or 7 < —1/2 (in which case the gas is localised). 
Thus there exists a critical interaction strength of 



and the single particle energy eigenvalues are 



2mR'^' 



k = 0,±1,±2 



(A2) 
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